Orbital current patterns in doped two-leg Cu-0 Hubbard ladders 
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In the weak coupling limit, we investigate two-leg ladders with a unit cell containing both Cu 
and O atoms, as a function of doping. For purely repulsive interactions, using bosonization, we find 
significant differences with the single orbital case: a completely massless quantum critical regime 
is obtained for a finite range of carrier concentration. In a broad region of the phase diagram the 
ground state consists of a pattern of orbital currents plus a density wave. NMR properties of the 
Cu and O nuclei are presented for the various phases. 



I. INTRODUCTION 

Over the past two decades, the description of strongly 
correlated electron materials has been one of the most 
actively pursued problems in condensed matter physics. 
When the strength of Coulomb interactions between car- 
riers is of the order of (or larger than) their kinetic 
energy, many new remarkable phenomena may occur. 
Their fingerprints are seen in experiments done on sys- 
tems such as cuprate compounds with high temperature 
superconductivity^ , cobaltites with large termopower-, 
magnesium oxides with colossal magnetoresistance^, or 
heavy fermions^. Among these materials, cuprates play 
a special role. At half filling they are insulators with an- 
tiferromagnetic (AF) order, but, with doping, a sequence 
of phases is observed including spin-glass, pseudogap, d- 
type superconductivity (SCd) and eventually Fermi-like 
behavior for very large carrier concentrations. 

Unfortunately there is, to date, no consensus on a 
theoretical model that would allow one to describe the 
physics of the Cu-0 planes. In order to get insight 
into this strong correlation problem, the study of lad- 
der structuresSiS, has proven quite useful. Ladders are 
the simplest systems that interpolate between one- and 
two- dimension. They constitute the quasi- one dimen- 
sional analog of the Cu-O sheets and, because of the 
reduced dimensionality, even weak interactions lead to 
dramatic effects. In the one dimensional (ID) case, the 
weak- and strong- interaction limits are usually smoothly 
connected^. Controlled non-perturbative methods - like 
bosonization or conformal field theory - and numerical 
techniques can be used to analyse these systems. 

Compounds characterized by a ladder structure^i^, 
such as SrCuO, have been synthesized. They show a 
variety of unusual properties, for example large magnetic 
fluctuations, SCd with purely repulsive interactions and 
metal-insulator transitions under high pressur o^'^'^°i^^'^^ . 
For these materials, increasing the pressure amounts to 
changing the bandwidth, and hence the ratio of Coulomb 
to kinetic energies in the ladder structure. 

These experimental developments provided a strong 
incentive for theorists to study two-leg ladders 
with Hubbard interactions between electrons. In 
the weakly interacting limit, renormalization group 



(RG) analysisi^ was used to explore their phase 
diagra m^'^d^d^d^, 18^19, 20 ^ Tsuchiizu and SuzumuraSi, 

and Tsuchiizu and Furusaki^^ performed an RG anal- 
ysis in bosonization language in order to explore the 
regime of dopings close to half filling. Using current al- 
gebra, where spin rotational symmetry was introduced 
a priori in order to derive RG equations, Balents and 
Fisher— i2£ established the phase diagram of two leg lad- 
der versus doping, showing that there was interesting 
physics at finite dopings. They identified a sequence of 
phases, labelled CnSm with n (m) gapless charge (spin) 
modes. Numerical DMRG calculations focused on the 
large U limi t^'^'^'^i^^ , the so called t-J approximation at 
half filling2&i21i2S. The relevance of interchain hoppings 
on the low energy physics was also adressed^Si^i^i^. 

The above mentioned papers all assume that, in the 
low energy limit, the Cu-O system can be reduced to 
an effective single orbital model. In the context of two 
dimensional (2D) cuprate materials, such reduction to a 
single orbital model was proposed by Zhang and Rice^ 
and it allowed one to derive phase diagrams for these 
systemsii^i. However this simplification was called into 
question, and it was pointed out that it is necessary to 
retain the full three band nature of the model in order to 
capture the important physics^^*^. This issue becomes 
particularly relevant when one examines the possible ex- 
istence of orbital current phases. Such phases were ini- 
tially proposed for the Hubbard modelSI. They were 
subsequently analyzed by various author o^'^^i'^^'^" , but in 
slave boson and in numerical calculations one finds that 
they arc unstable. For single band ladder models, con- 
trolled calculations appropriate to ID reveal that for spe- 
cial choices of interactions - which must include non lo- 
cal terms - staggered fiux patterns are stable. This phase 
breaks the translational symmetry of the lattice^i^. Ac- 
cording to some authorsii^, the 2D version of this state 
(the DDW phase) describes the pseudogap phase of the 
cuprates. An alternative type of orbital current pattern, 
which preserves the lattice translational symmetry, was 
advocated to describe the pseudogap phase^i^. It then 
requires using a three band model. Recent experimental 
data taken from neutron measurements^^' and polar Kerr 
effect^ would be consistent with the latter proposal, but 
more studies are clearly needed to fully corroborate this 
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scenario. 

Motivated by these considerations, Lee, Marston and 
Fjarestad^ (see also Ref. 113) generalized the system of 
RG equations written in current algebra language by Ba- 
lents and Fisher to study the Cu-O Hubbard ladder. 
Their work was however limited to the half-filled case, 
where umklapp terms dominate the physics, giving rise 
to Mott transitions. In a recent rapid communication**^ 
we outlined the method which allowed us to map out the 
full diagram of the Cu-O ladder as a function of doping. 

The aim of the present paper is to provide details of 
our derivation, and to present new results which are ex- 
perimentally testable. In our work, oxygen atoms are 
taken into account at each calculation step, which al- 
lows us to probe their influence. First, they lead to new 
types of phases compared with the single orbital case: 
a Luttinger liquid (LL) regime is found for a fine range 
of dopings and, in a broad region of the phase diagram, 
the ground state displays an orbital current pattern plus 
density wave quasi long-range order. Our study thus un- 
derscores the importance of including these additional 
degrees of freedom in the structure, in particular with 
regards to the existence and to the stability of currents 
patterns. Although our results have been derived for the 
specific case of ladders, they have potential relevance to 
the physics of 2D cuprate materials as well. Second, spec- 
troscopic tools measuring local properties, such as NMR, 
are predicted to give different signatures depending on 
wether they probe Cu or O sites. In the large U limit, for 
2D cuprates, it is believed that spin fluctuations on oxy- 
gen sites merely track those on the copper sites^. The 
advantage of revisiting the issue in a quasi-lD context is 
that one can monitor spin excitations on oxygen atoms 
both in the small and large U limits using bosonization 
techniques. This is done in the present paper for various 
dopings in the small U limit; we do find differences be- 
tween the NMR signal on the copper and oxygen atoms 
at low temperature, when gaps set in, but not at higher 
temperature in the Luttinger liquid (LL) regime. 

The paper is organized as follows: in Sec. |TT] we define 
the model including the interactions relevant to the low 
energy physics. In the continuum limit the quadratic 
part of the Hamiltonian is diagonal in a particular ba- 
sis, Bo- We give the relations between this basis, the 
bonding/antibonding basis -Bo/rr (relevant in the non- 
interacting case) and the total/transverse density basis 
_B_i (the most appropriate to write "backward" interac- 
tions). 

In section III we present a new method which allows 
one to set up the RG equations in the case of generic 
doping4^. One of its salient features is that it treats 

the rotation of Bq with respect to B^ during the flow. 

This effect needs to be taken into account in order to per- 
form all calculations properly. We list the resulting set of 
equations; their derivation is presented in the Appendix. 



The various flows and the resulting phase diagram are 
given in section IV. First we assess the impact of the 
additional degrees of freedom, hence we set all Coulomb 
interactions pertaining to the O atoms and direct intcrox- 
ygen hoppings to zero. Some of the results obtained in 
previous work^^ can now be checked using our improved 
RG method. In constrast with the single orbital case, we 
find an intermediate doping range where all spin and gap 
modes are massless (i.e a quantum critical line). Next, in- 
teractions involving oxygen atoms and hoppings between 
these atoms are introduced. We find that interoxygen 
hoppings promote a phase of orbital currents and we an- 
alyze its structure. Spin-rotational symmetry was not 
imposed a priori, but we checked that the required prop- 
erty was preserved during the flow. This provides a check 
on the consistency of our calculations. In the case of 
massive regimes the evolution of the gaps with doping is 
shown. We briefly examine the impact of umklapp terms 
which are present at half filling. 

In section V, we compute spin correlation functions, 
which allows us to derive the Knight shifts K and the re- 
laxation rates Ti for Cu and O nuclei. There are several 
improved features in our work. In Ref. f2ll spin-spin cor- 
relation functions were calculated in the low temperature 
limit, using Majorana pseudo-fermions for the spin part; 
this assumes that gap opening in the spin and in the 
charge modes occur at well separated T. In bosoniza- 
tion language, spin-spin correlation function are easily 
obtained in all cases. For instance, if one treats spin and 
charge density fluctuations on equal footing, one shows 
that the uniform part of the susceptibility approaches 
a quantum critical point as doping increases and that, 
at low temperatures in the gapped phase, the staggered 
part gives a different temperature dependence for each 
atom in the elementary cell. We also discuss physical 
implications of the orbital current phase. 



II. THE MODEL 

A. Hubbard Hamiltonian for Cu — O two-leg ladders 

We consider a two-leg ladder with a unit cell con- 
taining two Cu and five O atoms. Two edge oxygen 
sites are included because they would provide connec- 
tions with neighboring ladders (which are not considered 
in the present work). The hamiltonian of this system 
is divided in two parts: the kinetic energy of electrons 
moving on the lattice Ht and electron interactions Hint 

H = Ht + H,nt (1) 
The explicit form of the first, tight-binding, part is 
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^\.^^mja{b^n3a + &mj-l,<T) + H.C.] 



m^Cu 



) + H.c.]) 



(2) 



J. I 




FIG. 1: Energy contributions to the hamiltonian of the Cu-O 
Hubbard ladder; the figure shows two-unit cells. Subscripts 
in the a, b annihilation operators track the coordinate of the 
various atoms in a cell. 



Ht — ^( ^ <^Cunmjcr + ^ ^Onmja — 
j(7 m^Cu ni^O 

- Y t_L[aLj^(6m+l,i<T + 6m-lja)+H.C.])- 



where amjaibmja) annihilates holes with spin cr on a cop- 
per (oxygen) site, j labels cells along the chain and m 
labels the atoms within each cell, rimj^ — ciJnja'^rnja is 
the density of particles on site to, and we use here hole 
notation such that t, t±, tpp are all positive, e = eo~£Cu 
is the difference between the oxygen and copper on-site 
energies. 

LDA determined values^ of the parameters pertaining 
to SrCuO systems show that interladder hopping ampli- 
tudes are at least one order of magnitude smaller than 
their intraladder counterparts, so that the two-leg lad- 
der description is an excellent starting point for these 
compounds. Inside the elementary cell, t and tj_ are the 
dominant hoppings and their values are comparable. The 
difference between the electronic Cu- d and O- p state 
energies Ecu and Eq, is about 0.5t. There does not ap- 
pear to be ab initio determinations of Coulomb terms 
for SrCuO ladders, but from what is known for cuprates, 
we may estimate a local U of order 5< for the Cu sites, 
meaning a strongly interacting regime. In the follow- 
ing we will use constant values of the band parameters 
t ^ t± ^ 1 and e = 0.5, and treat the other observables 
{Ucu, tpp/t, Uo/Ucu, Vcu-o/Ucu) as tunable variables. 
In order to gain insight into the physics of the multiband 
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case, we analyze the above model using a renormalization 
group procedure in the interactions, i.e we assume that 
all of these are smaller than the kinetic energy; hence, 
the validity of the solutions cannot be ascertained, in 
the event when some of the interactions were to grow so 
large during the flow that they became of the order of 
the bandwidth. As was stated above, the experimental 
regime corresponds to a situation where Coulomb terms 
are sizable. Nevertheless the RG approach allows one to 
obtain a full analytical solution of this complicated prob- 
lem, and to make detailed comparisons with the physics 
of the one band system. Furthermore, for the case of 
the single band ladder one finds that the physical prop- 
erties in the weakly- and strongly- interacting limits are 
smoothly connected. We will come back to that point 
when we discuss our results. 



Eigenvalues and eigenvectors of the non-interacting 
part are simply obtained by Fourier transforming Ht- 
Since e is of order t, we neglect the non-bonding and an- 
tibonding higher energy bands which are mostly of p-type 
character, and this reduces the model to two lowest lying 
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bands crossing the Fermi energy. The Hamiltonian is: 

Ht = ^ea{k)nkaa (3) 

where a = 0, tt denotes the bands and a are spin indices. 
The operators corresponding to the eigenstates of Ht are 

Omfco- (or bmka) ^2^\na0.aka (4) 



ea{k) are the eigenvalues of Ht (the Cu-O distance is 
set to unity), and Xma are the amphtudes of the over- 
laps of the eigenvectors with the atomic wavefunctions 
in the unit cell. This defines the bonding (o) and an- 
tibonding (tt) eigenbasis B^j^. For ipp ^ 0, the o and 
TT energy bands are the two lowest, real solutions of the 
characteristic equation 



(e - eo(fc))(3ii + - eo{kf) - 2(1 + cos(fc))(-6ippi_Lt - i^(e - eo(fc)) + -itl^^t + eo(fc))) = 
(e - e^(k)){t\ +€^- e^{kf) - 2(1 + co&{k)){^2tppt - t^{e ~ e^(fc)) + tlp{e + e^(fc))) = (5) 

en = e 



Including tpp increases the values of the Xbia for the O 
atoms and makes the o and tt bands more asymmetric, 
but there are still only two bands crossing the Fermi en- 
ergy, so that the analysis remains valid. We note, how- 
ever, that the contribution of the oxygen p-orbital per- 



pendicular to the one participating in the Cu-O bonding 
increases as tpp grows larger, until, for tpp > 0.5t, it dom- 
inates that of the copper d-orbital. Hence, we confine the 
range of variation of tpp to 0-0. 5t. 

The interaction part in fermionic language is given by: 



Hint 



Uonmj]n.mj[ 4- 



m^Cu,n^O a .a' 



ja') 



(6) 



B. The continuum limit and bosonization 

We now express the Hamiltonian in bosonic represen- 
tation. The procedure is standard'^iSS, and we only outline 
the main steps here. We linearize the dispersion relation 
in the vicinity of the Fermi energy: 

^ J2 ^(l^Faa'^arqaO-arqa (7) 

|g|<Q raa 

r = ±1 denotes right and left movers, with momenta 
close to their respective ifci;-, Q is a momentum cut- 
ofi^. The boson phase fields denoted by 4>aa{x) are in- 
troduced for each fermion specie, aa contains spin and 
band indices, x is the spatial coordinate along the ladder. 
Fermionic operators are expressed in term of the bosonic 
field 4iaa{x) and 0„a{x) related to carriers fluctuations, 
by 

)) (8) 

where Ty^a are the Klein factors which satisfy the required 
anticommutation relations for fermions. These rja-a do 
not contain any spatial dependence and they commute 
with the Hamiltonian operator. They only influence the 
form of the order operator in bosonic language (through 



terms of the form rjcraVa-a') and the signs of the non-linear 
couplings through a F coefficient (the eigenvalue of the 
'na-aVcra'ila-a"'>lcra"' Operator). The operator is unitary so 
F^ = 1. This equality applies also to linear combinations 
of fields (change of basis); Following Ref. [2^, we choose 
F = +1 in the a+/— basis (see below). We also introduce 
the phase field 9„aix)] its spatial derivative Ilaaix) = 
dxQaaix) is canonically conjugated to 4>aa{x)- 

Now the hamiltonian may be rewritten using the above 
phase fields. The interaction term in the Hamiltonian 
can be split in two parts; one part only depends on the 
density of right and left movers, and gives - as does the 
kinetic energy - a contribution quadratic in the fields (pu 
and 9^, (where v labels the eigenmodes in the diagonal 
basis) of the form 

^0 = E/ ^[{uM{^iiuf + i^WA.?] (9) 

For the non-interacting system, one has K^, = 1 for all 
modes, and Hq is quadratic in the diagonal density basis 
which is simply (the momentum k^ associated with 
the rungs is either or tt). Another basis commonly 
used in the literature is the total/transverse one It 
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is related to Bott by: 



± 



V2 



(10) 



where /i stands for spin or charge depending on which 
density is considered. 

In general X in ([5]) is a matrix, the form of which 
depends on the basis in which the densities are expressed. 
For example if we use at the start of the calculation 
(the basis which diagonalizes the tight-binding part of 
the Hamiltonian) , we obtain 



ii-R- 







9oo 


9L ^ 


fL 




9L 




fob 


9L 


Vfo 


9L 


\ 9^ 


9-K7T 




Vf. ) 



(11) 



and 



ii-K = 



Vfo 
Vf. 



Vfo 
Vftt 



(12) 



where Vfo/t, 

,-L(ll) 



are the Fermi velocities in the o and tt bands 



and Qij are interactions between electron densities in 
the i and j bands, with perpendicular (parallel) spin. 



In order to express the hamiltonian in a Gaussian form 
(Eq. [5]), which is quite convenient for the RG calculation, 
we diagonalize K. This defines the Bo basis. In general, 
Bo is neither the bonding/antibonding basis -Bottj nor 

the total/transverse basis B^ We define the S matrix 

which describes the relative orientation of the Bo and 
B, bases: 



/ -Pi Qi 

"Qi Fi 

P2 Q2 

V -Q2 P2 



(13) 



One can express the parameters Pi and Qi with the help 
of angles a (for the spin part) and (3 (for the charge part): 



Pi = cos a + sin a 
Qi = cos a — sin a 
P2 = cos (3 + sin /3 
Q2 — cos /3 — sin /3 



(14) 



The remaining part of the interactions has a non-linear, 
cosine, form, in bosonization language. The most conve- 
nient basis to express this contribution is B^ and one 

findsii^ 



^mt(i) =~ficy drcos{2(t)s+)-cos{29c-)+giaJ drcos{2{(t>s+)-cos{2{9,-)~g2cJ dr cos(2(6lc-) • cos(2(0,_)-t- 
+94:a J dr cos{2(j)s-)-cos{29s-)+gi J dr cos{2(j)s+)-coa{2(j)s_)+g2 J dr sm{2(l}s-) sm{2(j)s+)+g\\c J dr cos{26c-)-cos{29s-) 



(15) 



where F coefficients determine the signs of the gi cou- 
plings (for instance, this gives minus signs for gic and 
<?2c)- We use the following notation: indices 1 to 4 re- 
fer to the standard g-ology processes for the left and 
right moving carriers, letters a to d correspond to similar 
processes, when the o and tt bands labels are used in- 
stead of the left or right labels. The relation between 
the gi couplings and the ones in ^ is given in Ap- 
pendix 1^ Note that in the quadratic piece, both 52 
and (74 -type terms need to be included, in order to prop- 
erly account for magnetic fluctuations^^. Examples of 
interaction processes are shown in Fig. [21 For instance, 
the two gid terms describe events where one right- and 
one left- moving fermion, both belonging to the same (0 
or tt) band, backscatter within that band. If we bosonize 
this contribution, we find two terms, gid cos(0i +4>2) and 
gid' cos((/)i — 4>2), instead of 51(2) ■ <?i and g2 correspond 
to the sum and to the difference of these "Id" -type pro- 




FIG. 2: Diagram showing some of the scattering processes. 
Blue (green) lines are for carriers in the o{tt) bands. This 
illustrates the notation used for the cosine-type terms 



and 



cesses respectively (51 = il±^ML^ — EidLzMl±"j^ 
(72 7^ when the O atoms are included. If the two bands 
were equivalent only the gi process would be present. 

In a standard Hubbard model, only spin perpendicular 
terms are present at bare level, and the last term in Eq. 
(fT5)) does not appear at the beginning of the flow. How- 
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ever, Bourbonnais pointed out^-S. that, during the flow 
towards the fixed point, additional scatterings involving 
electrons with parallel spins are generated by the RG pro- 
cedure. In our case, we are including a Vcu-O term so 
that, right from the start, our model contains interactions 
between carriers with parallel spin. Vcu-o gives rise to a 
non-linear cosine term while the other spin-parallel pro- 
cesses, which are generated by the RG procedure, give 
contributions to the various K . 

The g^a term has a non-zero conformal spin and gen- 
erates two extra couplings during the renormalization: 

H^nt(2) = -Gp J dr cos(40,_) - J dr cos(4^^,_) (16) 

These additional terms need to be taken into account, 
because they might become relevant when the other in- 
teractions scale to zero. 



III. THE RENORMALIZATION GROUP 
ANALYSIS 

A. Incommensurate filling 

We start from the quadratic part of the Hamiltonian, 
and treat the non-quadratic part (jlSp in perturbation, 
using a renormalization group procedure. We compute 



the corrections to the correlation functions to second or- 
der in g, and we incorporate them into the LL parame- 
ters K. However g terms are expressed in the B-\ basis 

while the quadratic part ^ is diagonal in the Bq basis, 
so the Pi{a,(3) and Qi{a,P) coefficients come into play. 
As a result, off-diagonal terms are generated in the K 
matrix during the RG iteration. At this stage. Bo is no 
longer the diagonal basis. In order to fix this, the Bo 
basis has to rotate during a renormalization cycle. In 
addition to the standard RG equations for the interac- 
tions, we need to find the RG flow of the angles a (for 
the spin density basis rotation) and (3 (for the charge 
density basis rotation). So, first we determine the cor- 
rections dKi,...,dK4, dBi2, dB34 that change the entries 
of the K matrix during the intial RG phase. Next, we 

go back to _B-| , using the transformation S~^; Since 

B^ is a fixed basis, the increments of the K matrix 

elements give the RG step corrections expressed in the 
_B-i basis. This new matrix is diagonalized by the oper- 
ator S{a -\- da, (3 + d/3) where the angles da, df3 depend 
on c?_B^_^+ and dK^_(^^_^_) (/z = c, s). The procedure is 
summarized in the diagram shown in Fig. [31 A detailed 
derivation is given in Appendix [Bl where, for the incom- 
mensurate case, we set all umklapp terms to zero in Eqs. 
(|B16p and we obtain the following set of differential equa- 
tions 
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( Ki II \ 

K-2 II 

K:, 

\ il Ki J 



RG 



I Ki + tlK, 



(2) 



\ 



II I) Ks + dKi ilBu 

\ II dBjj Ki + tlK, j 




1) K,- 6,.-,+ 



II \ 

B,,_rt4ilS,,_,t ffrt+(IA;+ / 



FIG. 3: Diagram showing the flow of the diagonal basis during 
the renormalization. The bottom row, shows the matrix in 

the fixed B-i ; the top row corresponds to the diagonal basis, 

used to write down the RG equations. 



dKi 
"df 
dK2 

dl 
dKs 

dl 
dKi 
dl 

dgic 

dl 
dgia 

dl 

dg2c 
dl 

dg\\c 
dl 

dgia 
dl 
dgi 
dl 
dg2 
dl 

dGp 

IT 

dGt 
dl 



Kl{Qlgla + Qlglc + P?Gl + P^gl + ^{gf + gl) + f{Pi){gig2))] 



[Qiigia + + Gi) - KiiP^gt + Pigi, + QiGl + Qfg^c + 2 + 5l) - /(A)(5i52))] 
Ql[9lc + 9lc + glc] 



<?ie • [2 - {PIK2 + PiK-^ + QjK, + QlK-')] - (.gi52c + 5ia5|lc) 
g^a ■ [2 - iPl{K2 + K-') + Ql{K, + R-'))] - g,,g\\, 
2 - {PiK-' + P^K, + QlK-^ + QIK2)] - gicgi 



52 



<?l,,(2 - (P^xr' + QlK2' + Pl^T' + QlKi')) - giagic 

gia{2 - liP^Ki + K^^) + Ql(K2 + K^'))) 

<?i • (2 - {K2 + Ki)) + PiQi{K2 - Ki)g2 - 75ic52c 

-52 • (2 - {K2 + Ki)) + PiQi{K2 - K^)gi 

Gp(l - {PlK^ + QIK2)) + gUPl{K, - K-^) + Q\{K2 - K-^)) 
Gt{l - iP^K^' + QlK^')) + gUPK'K, + Xfi) + Ql{-K2 + 

I 



(17) 
(18) 
(19) 
(20) 
(21) 
(22) 
(23) 
(24) 
(25) 
(26) 
(27) 
(28) 
(29) 



The equation giving the renormalization of g2 measures 
the influence of the O orbitals. The other two are conse- 
quences of the gia term. Note that P and Q depend on 
a and (3 (see Eq.([T3])), and hence they change during the 
flow. 



Additional renormalization equations for the rotation 
of Bo are 



dcot2a _ {{dKi - dK2)ta.\iAa + dBi2) 

dl ^ Ki - K2 

d cot 2/3 {{dK3~ dKi) tan AP + dBji) 



dl 



K3~Ki 



■ dl 
■dl 



-1 (30) 
(31) 
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where the equations for dBi2 and dB^^ are 



dB 



12 



dl 
CLB34 
dl 



PiQiiigl + al + Gf) - KiK^igl + gl + +gl + G^^j) - KM{P,)g^g2 



P2Q2{gl 



■glc 



r 



(32) 
(33) 



The function h is defined by 



0.25(Pf 



}\)y 



(34) 



Note that we introduced 51, (72 ~ the sum and difference 
of the gidS in both bands - because the renormahza- 
tion of g2 involves only gi and g2- The derivation of 
the renormalization equations in this case is presented in 
Appendix [B] 

As was found in Ref. 53 the interband scattering pro- 
cess (type c) renormalizes the Fermi velocities in both 
bands to a common value. The additional equation tak- 
ing this effect into account is: 



c?7 



= i{^-i){glc + 9lc + 9lc) 



(35) 



The initial value of this asymmetry parameter is 70 



W{a - where a - Vfo+v^ 



2-^ V" 2' ' "^^^^^ ^ 2Vfo 

Including this effect does not change our results, but 
it allows us to determine whether intra- or inter- band 



scatterings dominate for a given solution of the RG flow. 

If fully spin-isotropic interactions are present in the 
fermionic Hamiltonian, SU (2) spin-rotational symmetry 
has to be preserved during the RG flow. Some additional 
constrains on the RG variables can be derived in this case. 
For example one of them (for type c scattering) is 



g2c - 9ic - 9\\c = 



(36) 



Rather than using these constraints to reduce the number 
of RG equations, we check that they are satisfied during 
the flow. 



B. Half filling 

If the two-leg ladder is half-filled, additional umklapp 
terms should be included in the Hamiltonian 



Humk ^ g3\\ I drcos{24}s+)cos{2(l)c+ + Sx)+g3a I drcos{26s-)cos{2(j)c+ + 5x)+g3b J dr cos{2(l)s-) cos{2(l)c+ + Sx) 

+ gsc J dr cos{29c-) cos{2(f)c+ + Sx) (37) 
I 



Since these terms oscillate with S, their influence becomes 
important only for very small doping. The extended sys- 
tem of differential equations describing the RG flow has 
extra terms, compared with the incommensurate case, 
and each of them is multiplied by a doping dependent 
coefficient Jo{5). The full set of equations is given in 
Appendix [B] 

For small 6 these Bessel functions Jo{S) may be approx- 
imated by one and for large S by zero^i^. Starting from 
a small but non-zero doping, assuming that the chemical 
potential remains constant during the flow, the renormal- 
ization equation that describes this Mott physics is^ 



d5 
U 



{9' 



3|| 



+ gia + glb + gl)-MS) (38) 



The above equation gives an easy way to check if one is 
in the insulating or in the metallic phase, and which set of 
RG equations (with or without umklapp terms) is valid. 
S{1) flows to zero for the insulator and to infinity for 
the metal. The value of Sc depends on the initial values 
of gsi. The description of this transition is similar to 
that found in Ref. l3ll . which focused on the confinement- 
deconfinement transition of two-chain systems. 

For the sake of completeness, let us mention that other 
types of umklapp terms may appear for the two-leg lad- 
ders. These correspond to scattering of electrons in the 
bonding or antibonding bands, a process which becomes 
important if one of the kpi is around ^. In the presence of 
a large t± this condition may be fulfilled for dopings very 
different from zero. For tj^ > 0.1 1 it happens somewhere 
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in the C2S1 phase. As was pointed out in the discus- 
sion of the incommensurate case, couphngs involving 9c- 
flows then to zero. Thus both in the charge- and in the 
spin- sectors, one observes the rotation from the diago- 
nal basis to the k± — 0/tt basis. There are no processes 
competing with this, so the only effect is the appearance 
of a CIS"! region inside the C2S1 phase. These processes 
will not be considered in the following. 

IV. PHASE DIAGRAM 

Using the system of RG equations we determine the 
phase diagram. We identify the various phases based on 
the behavior of the renormalized quantities gi. We iterate 
the flow up to a point when some couplings become of or- 
der one. As usuali, the bosonized form is very convenient 
to analyze the strong coupling case, since when coeffi- 
cients in front of cosine-like terms become large, the cor- 
responding variables become locked. Subsequently, one 
may compute the physical observable in the ground state, 
by looking at the various order parameters in bosonic rep- 
resentation. These operators are given in Appendix ]C\ 
Some of the operators will now have exponentially de- 
creasing correlations, while others will decay as power 
laws. The dominant phase is the one for which correla- 
tions decrease with the smallest exponent. It corresponds 
to a quasi-long range order in the ladder. 

Two main factors may significantly affect the phase 
diagram that was predicted for two-leg Hubbard ladders 
with a single orbital per site: one is the asymmetry in the 
g terms due to the fact that the projections of the Cu 
and O orbitals onto the and vr bands have unequal am- 
plitudes and one is the influence of the extra parameters 
Uo, Vcu-o and tpp. 

We first investigate the impact of the asymmetry by 
setting Uo = Vcu-o — tpp — and we choose a small 
initial values for Ucu (in the range 10~^ 10~^). After 
this main part we consider a few additional issues such 
as the spin-rotational symmetry and the stability of the 
fixed points. 

As in Ref. we find that the parameter which de- 
scribes the behavior of the differential equations system 
is (3 = ^^2Vf^^'" • ratio ^ is constant, a depends 

only on 5: it is equal to one for half filling (then 5 = 0) 
and it reaches its maximal value when the Fermi energy 
is near the bottom of the band. The parameter a is only 
meaningful if the Fermi energy crosses both bonding and 
antibonding bands. We restrict our analysis to this case, 
otherwise one has a single band LL. 

A. Commensurate case 

Equations describing the commensurate situation are 
given in Appendix IB 31 

In this limit, umklapp terms lead to insulating phases 
with a gap in the charge degrees of freedom. These states 



are quite similar to those presented in Ref.[4^. In Sec. lVBI 
We will discuss this case and also similarities and differ- 
ences with previous studies. 

New and interesting physics occurs when the ladder is 
doped away from the commensurate case, and we focus 
on this situation in the remaining parts of this section. In 
the incommensurate case, the asymmetry that is present 
when the unit cell contains two different atoms (Cu and 
O) plays a critical role and leads to differences between 
the single- and multi- band models. 

B. The small doping case 

For small 5, cot 2a and cot 2/3 0, so the 
total/transverse density basis is the eigenbasis at the 
fixed point. In this case, 32, 54a, Gp, Gt are irrelevant. 
In the notation of Balents and Fisher^^, this is the CI SO 
phase where only the c-\- charge mode is massless. Fields 
Oc- and (j)s+ are ordered with the following values (given 
mod 2tt): 9c~ — 0, (f>s+ = 0. For s— the mode (spin- 
transverse), terms involving both 0s_ and 9s-, which 
are canonically conjugated become relevant, so one ob- 
serves an ordering competition between these two fields. 
The analysis of order operators presented in Appendix 
B shows that SCd-type fluctuations dominate if (j)s- is 
locked at 0, whereas if 9s- = an orbital antiferromag- 
netic state (OAF) is preferred. In our model, SCd always 
dominates for repulsive Ucu- This prediction confirms 
many previous discussions of SCd in two-leg ladders, in- 
cluding in the strong coupling regime^l, and in an inho- 
mogeneous doping situation^^. 

The advantage of working in bosonization language is 
that one can find a reasonable quasi-classical limit of the 
strong coupling fixed point. Using a semiclassical ap- 
proximation for the Sine-Gordon model^ allows one to 
find the doping dependence of the gaps in the system, 
which up to now was only obtained numerically . The 
following expression for the soliton mass (it is the lowest 
lying excitation if 0.5 < Ki < 2) is used 

where Ui, and K^, are the velocity and LL parameter of 
the v-th mode (by definition we are working in the di- 
agonal basis), and g is the interactions which makes this 
particular mode massive. For a more detailed discus- 
sion of gaps evaluation using RG see for example Ref. [H. 
The plot shows the behavior of the masses versus doping, 
evaluated with the above formula. One sees that, in the 
SCd phase, spin gaps go to zero as doping increases and 
so does the charge antisymmetric mode which has the 
largest value. The behavior of the gaps for small doping, 
showing a rather slow decay of their values is in agree- 
ment with experimental observation^i^. It is also com- 
parable with predictions obtained after refermionization 
of the problem and mapping it onto an exactly solvable 
Gross- Neveu model with SO (8) symmetry^ (but strict 
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FIG. 4: Doping dependence of the gaps in the spin transverse 
("1"), total ("2") and charge transverse ("3") density modes 
for the SCd phase. The rrii, are given in units of Ao, the initial 
energy cut-off of the RG procedure (of order ~ leV). 



duces only one massive spin mode, and we get the C2S1 
phase predicted by Balents and Fisher—; the interaction 
term which causes this behavior in our case is identical to 
theirs, once expressed in current density formahsm. Us- 
ing the same method as the one described for low dopings 
we are able to evaluate the doping dependence of the gap 
of strongly doped ladders. The result is shown in Fig. [S] 
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constraints for the Fermi velocities - viz Vpo = Vftt - and 
for the ratios of the g couplings at the fixed point have 
to be fulfilled then). There were also attempts to reduce 
the low-energy physics of pure Cu two-leg ladder to an 
SO (5) symmetric case^^iSi. For our more general system 
these conditions are usually not met, and spurious phases 
may even appear if one breaks some of the symmetries, so 
unfortunately we are not allowed to use these integrable 
models in our calculations. However some predictions, 
like the decrease of the gaps with doping and their rel- 
ative magnitudes are in complete agreement with these 
special cases. 

It is also worthwile pointing out that the values of the 
two spin gaps are always comparable, so that the ap- 
proximations that are made when mi >> m2 cannot not 
be used here to calculate the physical properties of our 
model. This behavior pertains to the range 6 < dd] upon 
approaching 6ci — 0.2 from below, gaps tend to close. As 
we will show next, a different phase, C2S1, emerges for 

5 > 6c2 ■ The intermediate range Sd < S < dc2 will be dis- 
cussed separately, when we examine the transition from 
the C15'0 to the C2S1 phase. 

C. The large doping case 

If the asymmetry between the bonding and antibond- 
ing bands is larger, cot2Q: —^ oo and cot 2/3 — oo, 
signalling that Bq-k is the eigenbasis for both the spin 
and charge modes. The RG flow converges very quickly 
to that fixed point for large dopings, typically when 

6 > 0.41. This corresponds to a = 4.2, a value that 
agrees with that found previously in Ref. [Tol . Interac- 
tions are not able to renormalize the ratio of the band 
Fermi velocities 7 to one anymore, which confirms that 
the BoTT basis is relevant for this regime. In the large 
doping phase, only gi ~ —172 are relevant. If one takes 
into account the rotation of the diagonal basis, which 
occurs as cot 2a varies, it appears that this flow pro- 
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FIG. 5: Doping dependence of the gap for the "o" band spin 
mode. Note that the gap is given in units of the cut-off Ao, 
the value of which decreases quickly when one approaches the 
bottom of the upper band 

One can easily identify the nature of the C2S'l phase, 
in bosonic fleld language: because the diverging interac- 
tion is go — gi — g2 <- 0, the slowest decay of correlations 
is observed for the CDW operator in the bonding ("o") 
band. 

For large enough dopings {S > 0.41), a gap will open, 
even if one starts from very small bare values of the in- 
teractions. In the range 0.28 < S < 0.41, angles still flow 
to the fixed point limits cot2Q;(/3) — s- oo, but gi and 
g2 grow very slowly. One needs to choose larger initial 
values of the bare interactions (but still smaller than the 
hopping t) and/or assume that one is close to the C2S'l 
region (starting from large cot 2a{P)) to find the gap ex- 
actly in the relevant spin mode. We conclude that the 
C2S1 phase exists in the entire range S > Sc2 = 0.28 and 
that, when 0.28 < 5 < 0.41, gi and 52 are very weakly 
relevant and thus very sensitive to higher order correc- 
tions. 

Previously, the existence of a C2S2 massless phase was 
predicted very close to the bottom of the band, (i.e when 
a becomes quite large). Our calculation, however, shows 
that the C2S'l phase remains stable in that limit. The 
reason for this difference stems from the choice of initial 
conditions in Ref. For single orbital ladders, when 
only on-site Hubbard interactions are included, the initial 
g2 is accidently zero. In our case, the presence of O 
orbitals always implies a non-zero initial g2- This says 
that the 32 term drives the transition, for very large a. 

At the bottom of the band the dispersion is quadratic, 
so the bosonization procedure, which requires a linear 
spectrum around the Fermi points, is not valid. The 
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calculation is done using conventional diagramatic tech- 
niques, and it confirms the stabihty of the C2SI phase 
with the same relevant coupling as found before. 



D. Quantum critical regime between the CISO and 
C2S1 phases 

We now turn to the intermediate regime 6ci < S < Sc2- 
Our key finding is that in this range, a new, massless, 
phase exists that had not reported previously, because it 
is found in the asymmetric limit, i.e when the unit cell 
contains both Cu and O orbitals. 

When one approaches the range S e (0.2; 0.28) either 
from below or from above, gaps appear to go to zero (see 
Figs. [3] and [5]). For Sd < S < Sc2, one has a line of crit- 
ical points where the phase is totally massless (C2S'2). 
Strong fluctuations, in particular near the critical end 
points 6ci and Sc2, cause poor convergence of the RG 
differential equation system. What is more, the angles 
Eq. (|30p vary significantly in a narrow range of I. Us- 
ing controlled approximations, we obtain an analytical 
solution that reveals the behavior of the system in this 
phase: we approach the critical end points from the mas- 
sive phases; we only keep the dominant couplings, which 
give us a simplified system of equations. Next we analyze 
the equations describing the angle rotations, and look 
for the range were derivatives become large, which takes 
place close to the fixed points. This gives us a condition 
for the divergence of cot2a{l3). Once the fixed point is 
known, one may simplify further the RG differential sys- 
tem. Now computing the RG exponent of each coupling 
is straightforward and enables us to find those couplings 
which remain relevant within the range of interest. 

Let us first consider 5c2 = 0.28. This point corresponds 
to the initial value cot2Q! =1. A numerical solution 
shows that the signs of {dKi — dK2) and B12 are the same 
and positive whereas the sign of {Ki — K2) is negative. 
From this simple analysis we infer that below this value 
I cot2Q;| decreases to zero and that above, it increases to 
infinity. Now, gi and (72 are only relevant when {K2 + 
Ki))—PiQi{K2 — Ki) < 2. needs to decrease strongly 
for this condition to be fulfilled and it is necessary to 
have nonzero values of f{Pi) and Qi at the fixed point. 
This condition corresponds to | cot 2a| 00, so one sees 
that below Sc2 = 0.28 the gi and 52 couplings cannot be 
relevant. 

The analysis pertaining to (5ci = 0.2 is less straightfor- 
ward. It involves | cot 2/3| and, because {dK^ — dKi) and 
_B34 have opposite signs, it is harder to get the fiow cor- 
rectly. The transition between going to zero and diverg- 
ing takes place when the absolute value of the two terms 
are equal. A detailed analysis of the angle dependent 
part of d cot 2/3 shows that this happens for 5ci — 0.2. 
When I cot 2/? I — > 00, then K^^ , which is much larger 
then one, influences the renormalization of the cos 6'c- 
coupling on equal footing with K^^ . This is the reason 
why these interactions are not relevant anymore. 



The above first order RG analysis was done in the 
vicinity of the critical end points and proves that when 
& E (0.2; 0.28) all interaction terms which are relevant 
outside are irrelevant inside this range. We have con- 
firmed the above simplified analysis by performing a 
numerical analysis of the full set of equations which 
shows that no other coupling is relevant. We see that 
a C2S2 phase is present between the CISQ and C2S'l 
phases. Hints for the possible existence of such state 
came from numerical studiesSS or from some special mod- 
els of ladders^^i^ with specific types of geometries, but 
we give here a direct proof of the existence of this phase 
for a generic ladder. 

The C2S2 phase is a LL where Bo-n is the fixed point 
eigenbasis for the charge modes and 5+ that for the 
spin modes. As far as the charge modes are concerned, 
K4=o is significantly smaller than one, while K3=t^ is very 
close to one at the fixed point. The spin parameters are 
both close to one because of the spin rotational symme- 
try. Thus one expects that correlation functions of band 
density fiuctuations of the form c^/_ ^criC-|_/_.o have the 
slowest decay. Logarithmic corrections need to be eval- 
uated, owing to the presence of a (single) marginal cou- 
pling go = (71 — 52 > 0. They show that a SDW within 
the o-band {SDW(o)) is dominant. 

E. The influence of Uo and Vcu-o 

Sofar, we have only discussed changes that stem from 
the presence of O orbitals in the structure. We now turn 
on the interactions involving the O atoms - Uo and/or 
Vcu-o ~ and probe whether these additional terms af- 
fect or not the phase diagram that we have found previ- 
ously. In the following, we assume that these interactions 
do not generate new types of terms in bosonization lan- 
guage, but that they modify the initial parameters of the 
fiow (for a detailed discussion of V-type terms, see for 
example Ref . [65l ) . 

In the CIS'O phase, SCd becomes less stable if large 
Uo or Vcu-o are present, but it has always a lower free 
energy than the OAF phase. When both Uo and Vcu-o 
are present, they seem to have competing effects. One 
would need to assume a very large attractive bare Vcu-o 
{Vcu-o < — 3.6J7c«) in order to stabilize a phase differ- 
ent from SCd. It would be a SCs phase with (ps- = 0, 
<f>s+ — 0, 6c- = 7r/2 and it would be very robust, even 
if the Fermi level approaches the bottom of the tt band. 
The existence of this phase, generated by Vqu-Oi was 
first pointed out in Ref. ,46i. The discussion of Ref. per- 
tains to the half filled case. The nature of the phase tran- 
sition between the two superconducting phases in ladders 
was described in detail in Ref. [2^ . so we are not going to 
discus this point. In the physical range of values of bare 
ycu-Oi one does not expect SCs to dominate. 

In the C2S'l phase, Uo and Vcu-o do not change the 
results significantly. Their main influence is that they 
make the gap smaller. It is to be expected, since the 
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CDW in the "o" band has an overlap with O atoms sit- 
ting between two Cu and introducing electron repulsion 
on O atoms makes the CDW less stable. For very large 
attractive Vcu-o the SCs phase re-enters. 

In the C2S2 phase, increasing Vcu-o has little effect 
on Sc2 but shifts 5ci to higher values. For ^'jj"^" > 1 
the quantum critical line still exists and an unphysically 
large ratio of « 5 would be required to suppress the mass- 
less phase and to observe a reentrant CI 5*0 phase with 
superconducting fluctuations. The Sci phase boundary is 
not affected by Uq or by Vcu-o- 

V. DISCUSSION AND CONSEQUENCES 

In this section, we discuss our findings in connection 
with previous work done on ladders. We also show that 
the C2S2 and C25'l phases possess an orbital current 
quasi long-range order and we compare our result with 
other proposals of current patterns for cuprates. 

A. Differences with the single band case 

In the derivation of the RG equations using current al- 
gebra, total particle-hole symmetry was assumed. Yet, it 
was shown^ that V-type interactions, for instance, can 
generate terms which break this symmetry at the begin- 
ning of the flow. They generate the following terms 

• a sine interaction term g2 ^ g^to ^ 9i-n 

• interactions such as — and g\g — g\^, of the 
form V(/)s+V0s_, which are included in the defini- 
tion of the non-diagonal part of the K matrix (see 
Eq. ([5])); this implies, for instance, that Pi ^ P2 

• ff4-type interactions generate different velocities for 
the spin- and charge modes {per se this is not a 
relevant perturbation but it enhances the impact 
of the other two contributions). 

When O atoms are included between Cu atoms, even 

A* 

if only Ucu is present, the bare go ~ y^Ucu is dif- 
x"- 

ferent from 5^ ~ y^Ucu and particle- hole symmetry 
does not hold anymore (in the limit E —^ 00, one has 
^ao ~ ^aTT ~ and Vp ~ so these two effects 

cancels out and go — g^ is still ~ 0{1)). This shows that 
it is quite important to include the oxygen atoms in the 
description of the two-leg ladder. 

The system of RG equations that we have derived does 
not impose such particle-hole symmetry constraint, and 
hence it may flow to a new fixed point which corresponds 
to the C25'2 phase. At the fixed point, is the diago- 
nal basis for the spin modes and Bott is the diagonal basis 
for the charge modes. It should be emphasized that for 
all other phases (which had been found previously for sin- 
gle orbital ladders), the diagonal basis at the fixed point 



is the same for the spin and for the charge modes. The 
presence of the three bands thus allows the symmetry 
between spin and charge bases to be relaxed during the 
fiow, and is instrumental in stabilizing the (725*2 phase. 
For the case of a single band, Ref. 31. pointed out that 
two additional considerations could lead to a significantly 
modification of the phase diagram obtained in Ref. [l9l . 
using a weak coupling perturbative approach. The first 
one was the inclusion of all interactions, not simply the 
relevant ones, the second one was the stability of the fixed 
points. For example in Ref. [s^, it was argued that the 
stability of the C2S'l phase was compromised, because of 
"a spin proximity effect" . However this C251 phase was 
found in DMRG numerical studies^. In our calculation 
it is important to note that all possible interactions were 
taken into account, and we did not impose any a priori 
symmetry. The presence of the (725*1 phase, that we do 
find in our calculation, is thus intimately connected with 
the rotation of the spin basis towards the fixed point Boj^ 
eigenbasis. 

At each step, we monitored the spin rotational invari- 
ance of the Hamiltonian through Eq. ((36l) to check that 
our equations were producing a reliable flow. The result, 
for the case of small as well as large dopings, is displayed 
in Fig. [HI In addition, in the (725*1 phase, there is one 
single massless spin mode, so that its K parameter must 
remain equal to one during the flow; we verified that this 
property does hold. 



B. Half filling and close to half filling 

At half filling, the charge symmetric mode becomes 
massive. All modes are gapped and spatial correlations 
decay exponentially. This is due to three relevant umk- 
lapp couplings. The spin and charge transverse modes 
are locked into the same minima as before, and the tran- 
sition only affects the total charge mode. One may view 
this transition as a quantum order-disorder Ising type. 
At half filling the dominant phase is the quantum disor- 
dered D-Mott phase, which, upon doping, turns into SCd, 
its dual counterpart. For large attractive V, an S-Mott 
phase, the dual counterpart of SCs, dominates at half 
filling. When we vary the strength of the interactions, 
the boundaries between these two phases look similar to 
those found for incommensurate fillings. 

The half-filled case for the Cu-0 ladder was discussed 
in Ref. [i^ both in the weak and in the strong coupling 
limits. For weak interactions, we may directly compare 
their results with ours. They used current algebra to 
treat the low energy physics of ladders with and without 
outer oxygens (five and seven atoms in the unit cell re- 
spectively). In the latter case, a spin-Peierls phase (BDW 
in Appendix. B) dominates, whereas a D-Mott phase is fa- 
vored, in the former case. The authors claim that this 
difference is due to a larger leg to rung anisotropy when 
outer oxygens are not present. The outer oxygens were 
taken into account in our model but we nevertheless find 
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FIG. 6: Evolution of g2c (yellow) and gic + g\\c (red) versus 
number of RG steps, during the flow. It shows that spin 
rotational symmetry is preserved (a) small a phase (b) large 
a phase 



in the literature. 

For instance in the case of Cu-O ladders close to half 
filling (the strong coupling case discussion in Ref. a 
t-J approximation was used. It gave a uniform phase - 
related to D-Mott - in a broad region of positive Ucu- 
Vcu-o phase space. For large attractive Vcu-O, a phase 
with holes localized around copper atoms is found, prob- 
ably connected to our SCs ordering. Our phase diagram 
matches the above description. The SCd phase, which 
we find close to half filling, is clearly seen in numerical 
studies. The CISQ phase was connected with this type 
of ordering in Ref. [23 were it was also shown that the gap 
preventing a DW-type ordering decreases upon increas- 
ing the doping. In Ref. it was found that the region 
where this phase is stable can be extended up to U=4t. 
In the t-J model, the rigidity of SCd with respect to a 
finite difference in the chemical potential of the two-legs 
was also established^. The same type of ordering (rung 
singlet) also dominates at half filling, for a special choice 
of parameters giving an SO (5) symmetry2£, since in that 
case one may solve the model exactly. All these results 
were obtained for single band ladders; quantitative differ- 
ences occur when oxygen atoms are included in the unit 
cell, and these were analysed in a numerical study^. 

Few studies were devoted to the intermediate and 
large doping regimes; we discussed the C2S1 phase (see 
above), and, as far as the C2S2 phase is concerned, a 
DMRG study^^ suggested the existence of such gapless 
phase well inside the bands for a zig-zag ladder; the oc- 
currence of a massless phase in the strong interaction 
limit would be certainly remarkable. 



Orbital current patterns at intermediate and 
large dopings 



a D-Mott phase. More generally, our entire phase dia- 
gram is very similar to their "five orbital" case. A possi- 
ble reason for this discrepancy could be that 

their tpp is barely less than the Cu-O hopping ampli- 
tudes. In our calculations tpp is much smaller, in accor- 
dance with LDA studies^ and with experiments. Recall 
that, as was described in Section II, whenever tpp > 0.5i 
(or t_\_) non bonding p-orbitals become relevant degrees 
of freedom, and these were not included in our model. 
Similarly a large initial value of the nearest-neighbor in- 
teraction V causes an exchange of the weight of the d and 
p orbitals in the lowest lying bands during the RG flow. 
This limit is beyond the range of validity of a simple RG 
approach. 

The strong coupling case (the so called charge-transfer 
regime) is important, because, for real inorganic mate- 
rials, U is usually of order 5t. Still, two features of the 
weak coupling regime remain valid in strong coupling: 
one is that spin-charge separation holds and two is that 
in the SCd phase (for instance) there is still an exponen- 
tial decay of DW operators. A connection between the 
phase diagrams of these two regimes is often suggested 



In the previous sections, we had set tpp = 0. We now 
assess the influence of this hopping term on the phase 
diagram. As long as tpp < 0.5t, our RG method remains 
valid, and only the initial parameters are changing with 
tpp. As we increase tpp, we note that a increases for a 
given doping, but that both Sci and Sc2 are decreasing. 
For tpp — 0.5, their values are about half that quoted 
for tpp — 0. A negative tpp (the electron-doped case) 
has the opposite effect. This influence of tpp > can be 
understood as an increased asymmetry between bands. 

The phase diagram, which summarizes our study of 
the Cu-O ladder for carrier concentrations between half 
filling down to the bottom of the "upper" band, is shown 
in Fig. El 

A spectacular effect of tpp is that it leads to new types 
of current loops, involving oxygen sites; the range of pa- 
rameters where orbital current patterns (OCP) dominate 
is seen in Fig. [T] 

A finite tpp allows direct current flows between oxygen 
atoms, giving rise to additional patterns, enclosed inside 
the elementary cell. One of these preserves the mirror 
symmetry on the a axis (the axis parallel to the chain di- 
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FIG. 7: Phase diagram of two-leg Hubbard ladders versus 
doping for Ucu > 0. 5 — corresponds to the half filled case; 
umklapp terms which open up a gap in the charge symmetric 
mode are not included here; DW+OCP denotes phases with 
orbital current pattern (OCP) on top of a spin- (SDW) or a 
charge- (CDW) density wave. 
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FIG. 8: Current pattern in the C2S2 phase. It has a mirror 
symmetry with respect to the a axis plus an additional SDW. 
The modulation has an incommensurate spatial periodicity 



rection and passing through the mid- rung oxygens), and 
it is of special importance. This is because we have shown 
that, at least for moderate dopings, operators in the o 
band dominate. The current operator between two atoms 
"a" and "b" is defined as ja-h = Y.a,a' V'L^a'b-V'l'bV'aa 
(we sum over band indices) and the total current pattern 
operator is given as a sum of currents on each bond. For 
symmetry reasons, if the current pattern has a mirror 
symmetry along cr, then the total current operator has 
the form of a particle-hole fluctuation in the o band. 

Two conditions must be met in order to get a dominant 
contribution: the pattern must form closed loops origi- 
nating from and ending at Cu atoms and it has to possess 
a mirror symmetry with respect to the plane containing 
a and perpendicular to the plane of the ladder. One of 
the current patterns preserves both of the required con- 
strains and, since it is similar to a configuration proposed 
by Varma, we call it " Varmal" -type (Fig. [5]). Then, be- 
cause the current operator has the same dependence on 
phase fields as the DW operator, these fluctuations have 
the same power law decay. Computing their amplitude 
will tell us which type of order dominates. 

In the large doping regime ((7251), we compare the 
ampfitude of the "normal" CDW and of the OCP-fCDW. 
The amplitude of the latter is found by summing cur- 
rent operator contributions for loops with one Cu and 
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FIG. 9: Current amplitudes within the o-band, given in tpp 
units. For this particular doping, the pattern has an 8-cell 
periodicity; The figure shows half a period (the other half 
is simply obtained by repeating the amplitude pattern and 
reversing all the signs) 



two O atoms. We use the mirror symmetry and add 
first equivalent pairs of currents. Each of these pairs 
gives a contribution proportional to tijlm{\aa\ia) or 
tijlm(XiijaXbia), where tij is the hopping parameter be- 
tween the relevant atoms. We emphasize once again that, 
in the single (Cu) orbital case, Xaa = 1 so that the 
current operator between Cu atoms has the usual in- 
terband form. It is the presence of oxygens that gives 
Im(XaaXbia 7^ 0, allowing the geometry of a Varma-type 
pattern to appear in the theory. 

A numerical calculation shows that these quantities 
are of order one and change only by a few percent when 
the doping increases from 0.25 to the value of S at the 
bottom of the band. The result of this procedure (the 
amplitudes of the currents determined by the products 
of Xij coefficients) is shown in Fig. |9l since it is easier to 
visualize a commensurate pattern, we chose 6 = 0.9 such 
that kpo = 1/4 (note that only the "0" band would cross 
the Fermi level for such a large value of the doping). 
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Due to current conservation, the weakest link between 
atoms determines the maximum value of the current. It 
is clear that the magnitude of tpp determines whether 
or not the OCP+CDW state may exist. Since the total 
amplitude is proportional to 2tppJ2ij InT-i^bja^Ma) mul- 
tiplied by the number of links, it is straightforward to 
obtain a threshold value tp™*"'' « O.St above which the 
OCP+CDW phase dominates. 

Varma's work was concerned with the strong coupling 
regime in 2D, and the stability of the current patterns was 
studied in mean-field theory. The fact that we were able 
to find such a state in ID, in the weak coupling limit and 
with purely repulsive interactions, gives an interesting 
perspective on the possible existence of such orbital cur- 
rents. Note that a type ordering similar to the one we find 
(current 4-DW) has been suggested in numerical studies 
of two-leg ladders^ We will return to this issue in the 
last part of this section, where we make a contact with 
strong coupling results. The statement about the exis- 
tence of OCP states given above, obviously holds also for 
the (725*2 phase, where CDW fluctuations are replaced by 
SDW fluctuations (Fig. [8|). We also note two differences 
between our orbital current states and Varma's: in our 
case, the structure is incommensurate (the modulation is 
doping dependent) and we get an additional DW modula- 
tion. Hence, the OCP-I-DW state also share similarities 
with the DDW phased. The main point is that intro- 
ducing tpp gives the possibility of new types of current 
loops involving oxygen sites. Phases with time reversal 
symmetry breaking have been widely investigated, but, 
for single orbital models, currents flow along Cu square 
plaquettes, giving rise to the OAF state. As was shown 
in detail by Fjaerestad and Marston^^ they are described 
by inter-band creation-annihilation process ~ iplip-rr ■ The 
order operator in bosonization language is given in Ap- 
pendix D (it is the OoApif) operator) and, this type of 
quasi-long range order is stable provided one introduces 
an attractive Vcu-Cu- 

What about current patterns in the strong coupling 
regime? This issue was investigated numerically. One 
paper^^S showed that if time reversal symmetry is artifi- 
cially broken by adding a magnetic field, one promotes a 
state with OAF currents and a CDW modulation for the 
two-leg ladder, very similar to the one suggested in Sec. 
IV. Another— considered variants of t-J models, in hopes 
of finding current pattern phases. Although somewhat 
artificial values were assigned to some of the parameters, 
this study suggested that quasi long-range order of the 
current patterns could be obtained provided one changed 
the internal description of the rung. Furthermore, the 
current pattern is accompanied by a charge density wave 
structure. 

It should be pointed out that both papers established a 
direct connection between the strong and the weak cou- 
pling regimes. For example, Ref. \^ showed that the 
spatial decay of current-current correlations is similar in 
regions of parameter space corresponding to weak or large 
coupling RG. This paper also emphasizes that, in order 



to obtain current patterns, one needs to go beyond the- 
ories using properties of SO (5) symmetric models. This 
confirms our findings, that new physics emerges in for- 
malisms where symmetry breaking is a priori allowed. 



D. Experimental systems 

Experimentally, it is rather difficult to vary the doping 
in ladder compounds, and the large doping regime is still 
inaccessible^^. Furthermore, different methods (NMR, 
optical conductivity and X-ray measurements) yield dif- 
ferent values of the doping for a given systeroS. One 
of the most interesting compounds, Srii-xCaxCu2iOii, 
contains both chains and ladders, and it was shown that 
a change in pressure may cause a charge transfer between 
the two?. Calcium is also a factor that affects the car- 
rier content of the ladder. In the low doping regime, 
this system displays spin-gaps, as is well established in 
many NMR studies; this will be discussed in detail in 
the next section. As far as charge degree of freedoms are 
concerned, the situation is more complicated. There are 
optical conductivity measurements showing CDW order- 
ing in these systems at ambient pressure^. This kind of 
ordering may be due to a large Vcu-Cv^ not taken into 
account in our model or to inter-ladder electrostatic in- 
teractions. The SCd phase appears under pressure, with 
a maximum temperature of order lOK for an optimal 
pressure of 3.5 GPa. The role of pressure in this tran- 
sition is not clear: it may change the bandwidth, the 
couplings between ladders, the screening of the intra- or 
inter- ladder interactions, or the doping. Recently2&i21, 
soft X-ray measurements were performed for this sys- 
tem. Their main conclusion is that an insulating " hole 
crystal" phase exists for commensurate fillings. It is sug- 
gested that this phase melts for other dopings. The au- 
thors interpret their findings by invoking strong on-rung 
hole pairing. This analysis supports the picture that 
emerges from our study of the low doping regime. 



VI. NMR PROPERTIES 
A. Spin susceptibility and NMR relaxation rate 

The spin operator with momentum q is defined as 

^IT. ^"-1 + q)'^'c„.2(k) (40) 

k(Tlo-2 

where c = a,b (respectively the annihilation operator of a 

hole on Cu or on O) and cr* is a Pauli matrix. From linear 
response theory, the time-ordered susceptibility reads 

Xm'm(q: ^^n) = ^ / {Tr S'„,, {q, r)5;„(-q, 0)) exp(iw„T) 

(41) 
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The above function is defined only for Matsubara fre- 
quencies LJ„; taking the analytical continuation, one 
obtains the retarded susceptibility Xm'm~ ^^'^ hence 
derive^ analytical expressions for the measured NMR 
properties of the system. 

The NMR signal comes from a contact interaction be- 
tween a nucleus and the surrounding cloud of electrons 
in an s-orbital state. 

The temperature dependence of the shift in (Zeeman) 
frequency of the m-th nucleus stems from hoppings of car- 
riers from the m-th atom s-orbital to the highest occupied 
molecular orbital p ot d orbital of the neighbouring sites. 
Thus, the Knight shift is 

^™ = E ^'"'^(P - 0, c = 0) (42) 

m' 

where the summation is taken over all neighboring d-Cu 
and p-0 orbitals. The overlap coefficients Xma, which 
enter Xm,m', are evaluated using first order perturbation 
theory; we include hoppings between a s-Cu orbital and 
a, p-0 orbital on the neighbor sites or a d-Cu orbital on 
next-nearest neighbor sites, as well as hoppings between 
a s-0 orbital and a d-Cu orbital on neighboring sites. 

The spin-lattice relaxation rate is also affected by the 
electronic environment. The signal measured on the m-th 
nucleus is given by 

Tim Imleh^P ^ l^Zm 

In the following, we omit the "i" subscripts, because 
we are working with spin-rotationally invariant models. 
Taking into account the fact that the Fermi surface con- 
sists of pairs of points of the form i/cj?, the sum in ^ 
can be divided into two independent parts: a uniform 
piece (q around = 0) and a staggered piece (q around 
= 2kp). 

Using the Xma allows us to connect the time-ordered 
correlation functions Ra{r{x,T)) of carriers in band a 
(they are introduced in the bosonic phase field language), 
with the i?„i',„(r(a;, r)) defined for a site basis 

Rm'mir{x,T)) \Xm'o\'^\Kio\'^Ro{r{x,T))-\- 

|2i?,(r(x,T)) (44) 

In order to get the retarded ^,(q,uj) entering Eqs. 
(|42|43p . we use the fact that correlations for spin opera- 
tors and for their complex conjugates are equal, and we 
simply obtain the retarded spin susceptibility by a Wick 
rotation^^: 

X^^„(x, t) = 2e{t)Im[Rm'm{r{x, T))]r^^t+s (45) 

followed by Fourier transforming the last function. 

Because of conformal symmetry in our ID quantum 
theory, results for zero temperature correlations can be 



extended to finite temperatures by simply substituting 
for the complex coordinates the following expression 



r^[x,T,(3) ^ — lsin\Y[ ) smh( — ^-^ — ) (46) 

y 7r TT 

This substitution r^(a;, r) r^{x, t, /3) gives us the tem- 
perature dependence of the susceptibilities. 

This procedure is valid both for the uniform and for the 
staggered parts of the magnetization. We write the time 
ordered correlation functions Ra{r{x,Ty) in each band 
in terms of diagonal modes R,^{r{x,T)) — F[(/)^] for the 
staggered and the uniform part, separately. The form of 
depends on whether the v — th LL mode is massless 
or massive and it will be presented below. Given 
the substitution [46] allows us to obtain the temperature 
dependence of and Tj^, but as the temperature in- 
creases, the form of F[(j)^] is changing. Generally it is 
assumed that above the temperature corresponding 
to the value of the gap A^, thermal fiuctuations make 
the u-th mode massless. For example, at T = 0, in the 
CISQ phase, one starts with three gapped modes (two 
for the spin and, one for the charge); we increase T until 
the energy of the first gap As+ is reached. Above the cor- 
responding temperature Ts+ we may consider that there 
is effectively one gapped and one gapless spin mode, and 
similarly for the charge sector. The others gaps (A^- and 
Ac-) will successively close at temperatures Ts- and T^-. 

B. Doping dependence of the NMR signals 

A number of paper a^°'^^i^^ ,were devoted to the compu- 
tation of magnetic properties of two-leg ladders assum- 
ing symmetry entanglement at the fixed point {S0(5) or 
SO (8)). Yet, following the discussion in section IV, we 
use simpler, approximate, methods which nevertheless 
have a wider range of validity. 

1. Uniform part 

For the uniform magnetization, only spin correlations 
need to be taken into account. Because the spin den- 
sity is generally related to the spin phase field u{x) — 
■^V (j)a{x), the zero momentum part of Ro is a linear 
combination of bosonic correlations calculated in the di- 
agonal basis Ri, (r) . In the massless case it is known from 
LL properties and given by 

R. = ^ (47) 

i 

The contribution to NMR of these power laws has been 
evaluated many times before in the literature. One gets 
a dependence for the Knight shift and for the re- 
laxation rate. One may improve these result both in the 
high- and low energy limits. At low energies, logarithmic 
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corrections from relevant and marginal couplings (?(A) 
(where the energy scale A may be related to the temper- 
ature) should be taken into account. Then, gic and g2c 
contribute to the o and tt bands, gi + g2 to the o band 
and (?i+.g2 to the tt band. At high energies, the curvature 
of the bands may be taken into account using an RPA 
approximation, following Refs.— i^. 

In the massive case we use the massive Gaussian model 
to obtain fluctuations around the static quasi-classical 
solution (equilibrium position) and this leads to 

Ry = Ki{KQ{miri) + K2{miri)) (48) 

where we have used the fact that for harmonic fluctua- 
tions around the soliton 54>i{x) = 4>i{x) — (j)oi{x), corre- 
lations are given by a Bessel function {S(l)i{r)6(j>i{0)) ~ 
K(){mir). A first order expansion, valid for large r. 



gives exactly the same expression as that found in ex- 
act calculations^iiS^. One needs to evaluate the follow- 
ing integrals (the exact formulas for the LL x{k,uj) are 
known^^ but it is not necessary to use them here) 

J— ^ f dty2F[X^fixo)]R{xo,t) (49) 

J-lm J ~^ 



K,n = A,„, (fc = 0)(^A™, (fc = 0)) J dtdxR{x,t) 

in' 

(50) 

where the summation over xq accounts for the momen- 
tum dependence of \mv{k). For the uniform part, inte- 
grals can be calculated analytically 
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(51) 



the appropriate bound b,c is chosen for the Knight shift 
or for the relaxation rate, and depends on whether one 
integrates over a time or a space-time domain. 

The results for Knight shifts, calculated for different 
atoms and different dopings, are shown in Figs. [10] and 
[TT] The discussion of the relaxation rates is postponed 
until after the evaluation of the staggered part, because 
the quantity which is measured in experiments is the sum 
of the uniform and the staggered parts of the relaxation 
rate. 

For the (7150 phase, an activated behavior exp(— ■^) 
is seen for the Knight shifts of all the atoms. This shows 
clearly on the logarithmic plots shown in the inset of 
Fig. [TU] However let us stress that in the C150 case we 
have two spin gaps, so one expects a more complicated 
shape than a simple straight line. For higher tempera- 
tures the Knight shift saturates to a constant value. As 
is expected for the uniform susceptibility, the responses 
of the different atoms are similar, only their amplitudes 
are different (this is because of the |Am!/P coefficients). 
For larger dopings there are less electrons in the conduc- 
tion band, and their velocity is smaller, so the saturation 
value also decreases. As the doping increases, spin gaps 
decrease and curves saturate at a lower T until we reach 
the quantum critical point (QCP) at 5 = Sci- This be- 
havior for the susceptibility was described in Ref. [2l|. In 
their case, a QCP appears in the presence of Vcu-Cu] in 
our case, doping drives the transition. 

For the (725*1 phase, we obtain a finite susceptibility 
even at zero temperature. This feature comes from the 
massless spin mode. The central oxygen atom which is 



only coupled to the gapped band does not give a finite 
susceptibility at T=0. For the second, massive, mode we 
observe a behavior similar to the one described above for 
the CI 50 phase, with the single activation gap shown on 
the inset of Fig. [TUb for two dopings. 

For the intermediate doping C252 phase, including 
logarithmic corrections is the only way to generate some 
weak T dependence. They arise mainly from the pres- 
ence of the marginal gi and g2 terms. Their influence 
on the uniform susceptibility was described in detail in 
Refs<^i22. Differences in the amplitudes of the Knight 
shifts for the various atoms in the elementary cell stay 
pretty much the same from one phase to the next, since 
these amplitudes are simply determined by A coefficients. 



2. Staggered part 

For q = 2fci?, both the spin and charge parts con- 
tribute to the band correlation functions. The band -Rq/tt 
with 2kp wave vector is a product of a spin and and a 
charge part, Ra = Ra ■ Ra- The form of depends 
on the fixed point eigenbasis for the angles and on the 
possible existence of gaps. 

The expression for the gapped spin phase was ob- 
tained using the expression for the 2k p part of the 
spin density operator correlations which is given by 
{Oo/nSDwir)Oo/^SDw{0)f'"' cos(0i ± (/-a)- The last 
form could be evaluated using the fact that = 4'oi + S4>i 
where the fluctuations of S(j)i are described by a massive 
Gaussian model, as was shown in the case of the uniform 
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(a) 



FIG. 10: Temperature dependence of the Knight shifts for (a) 
the different atoms in the elementary cell and (b) different 
dopings in the (715*0 phase. The activation gap at low T is 
shown in the inset 
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For gapped spin (cr) modes there are two possibili- 



^ smh{K2Koim2r2)±KiKo{miri)){m2a)^^{mia) 
for the CISO phase 
R'^ = smh{K2Ko{m2r2)){m2a)^^ 

for the C2S1 phase (52) 

In the gapless case, one gets a power law behavior; 

for the high-T limit of the CISO phase, where the 

eigenbasis is relevant, we find 



Ra 



(53) 



where a corresponds to the band index and i = 1,3, 
j = 2, 4 are the LL modes; 

For the CISQ phase, the charge mode is only partially 
gapped: the field 9^ is locked so the charge antisymmetric 
mode does not give any contribution to SDW, but the 
massless , "4" (charge symmetric) mode gives a power law 
contribution 
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For the other phases, both charge modes are massless and 
in this case, -Bqtt is the fixed point basis, and we have 
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FIG. 11: Temperature dependence of the Knight shifts for the 
different atoms in the elementary cell (a) in the C2S2 phase 
and (b) in the (725*1 phase 
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(56) 



For the spin part in the (7251 phase, one substitutes 
K\ , K2 to Kj, , ^4 . The amplitudes of Ra on different 
atoms need to be calculated. Once again Xma{k) are 
involved, however for those atoms with neighbors along 
the ladder (on-leg Cu and O atoms) these coefficients 
are different, because of phase factors at k — 2k p which 
cause cancellations in some contributions of neighboring 
atoms. Another possible factor may cause differences be- 
tween atoms in the elementary cell. Following Ref. [H, 
one may assume that, below a characteristic distance 
X < Lg — , umklapp terms are relevant and that they 
open up a gap in the charge symmetric channel. This 
massive charge correlation affects the staggered part of 
the magnetic susceptibility, and yields an expression sim- 
ilar to Eq. ([5^ (with cosh instead of sinh). For on-leg 
oxygens, which sit between two Cu along the ladder, one 
recovers a sinh instead of a cosh. This produces different 
(54) amplitudes i?^,^(i) for Cu and for O on-leg atoms, pro- 



vided Lg > TO4 ^. We have made the calculation for the 
half filled case, and the result is that 777,4 is of the same 
order as 7772, so for dopings larger than 0.05 this effect 
should not play any role. 

Once band correlation functions are known one may 
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follow exactly the same procedure as in the uniform case 
in order to obtain the temperature dependence of t^-. 



3. Total relaxation rate 

The plots in Fig . [12] shows 1 /Ti for different atoms 
in the elementary cell. They were obtained by numerical 
integration of Eq. (|52p and adding the result to that 
computed for the uniform part. 

As for the Knight shifts, the difference between atoms 
are caused mainly by the different Xma coefficients. How- 
ever these coefficients can be different for the staggered 
part and for the uniform part. 

The first observation is the linear dependence of ^ 
at high temperatures for all atoms, for all dopings. The 
Knight shift saturates in this temperature range to a con- 
stant value and this is in accordance with the Korringa 
law. For the C2S2 phase we observe the linear depen- 
dence as expected for a massless LL with all K parame- 
ters close to one. 

The second main conclusion is that processes involv- 
ing large fcy transfers can strongly affect the measured 
rates, especially for temperatures comparable with the 
spin gaps, as previously reported^^. One observes only 
small differences between atoms in the elementary cell 
at low T. The difference in the relaxation rate of a Cu 
nucleus compared with a central O nucleus comes form 
the fact that the latter may only relax through processes 
in the "tt" band, while, for the former, both bands con- 
tribute. The difference between on-leg atoms and atoms 
sitting at other locations comes from the fact that the 
staggered part contribution of the former nucleus is very 
small, as it is suppressed by the opposite contributions 
of the two neighboring Cu atoms. The low-T activation 
behavior (in the C2S'l and (715*0 phases) is then clearly 
seen on these on-leg O sites. 

Two points should be kept in mind when comparing 
our results with experiments. First our Aq is of order 0.5 
e V, so that the largest charge antisymmetric gap is of or- 
der 700K; observing it would be experimentally challeng- 
ing, and it would be even harder to reach the Korringa 
regime predicted at higher T. Second, our calculations 
were made in the phase where SCd fluctuations domi- 
nate. Thus experiments done at large pressures would 
be the most relevant to compare our findings with. 

Aside from the above caveats, the results of our cal- 
culations seem to be in very reasonable agreement with 
experiments^ii^. 



VII. CONCLUSIONS 
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FIG. 12: Temperature dependence of the relaxation rates for 
the different atoms in the elementary cell (a) in the CISO 
phase, (b) in the (725*2 phase, (c) in the (725*1 phase 



differences between models which include O atoms and 
models which do not, even close to half-filling. The mass- 
less (7252 phase is of special importance in that respect. 
A Varma-like phase with incommensurate orbital cur- 
rent patterns and additional density wave characterize 
the ground state structure at intermediate and large dop- 
ings. Signatures of these states can be seen in NMR ex- 
periments probing the various nuclei in the cell. 

We see important differences between Cu and Cu-0 
ladders in the weak interaction limit (U < e), but nu- 
merical approaches which can investigate the opposite 
limit as well suggest that these differences do survive for 
U > e. This invites further analytical studies of the two- 
leg Cu-O Hubbard ladders in the large U limit. 



Our study has clearly shown that including oxygen 
atoms in the structure produces significant changes in 
the ground state phase diagram of doped, Cu-O, two- 
leg Hubbard ladders. This result is fully consistent with 
DMRG studies suggesting that there are quantitative 
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APPENDIX A: COUPLING CONSTANTS 

Initial conditions for the non-linear terms are 

= 4/(1, 2, 1,2,1) 
ff?, = 4/(l,2,2,l,l) 
52°e = 4/(1, 2, 1,2,-1) 
5f, = 4(/(l,2,l,2,l)-/(l,2,l,2,-l)) 

2/(1,1,1,1,1) ^ 2/(2,2,2,2,1) (Al) 



9i 



Vfo 



„^ 2/(1,1,1,1,1) 2/(2,2,2,2,1) 

54°a=4/(l,2,2,l,-l) 

where the function f{k,l,m,n,p) converts the interac- 
tions given in the atomic basis {Ucu,Uo,Vcu-o) into 
band "g-ology" interactions. 

fik,l,m,n,p) - ^ A,:,.A*;A,„,A,„(V^^"*™+ 

F™*^"- COs(fcj.,„ - p ■ kpn)) (A2) 

The summation is taken over all the atoms in the elemen- 
tary cell, y™*'"" denotes interactions within the elemen- 
tary cell, and is Vcm-o, since one of the atoms is 
outside the elementary cell, as in Ref. |4^. Initial values 
for K and cot(Q!), cot(/3) are evaluated as follows: start- 
ing from Eq. p^ . one performs a S{ti/A) rotation. In 

this basis, K can be calculated by simply solving a 

matrix equation. The initial K^, are given by the eigen- 
values of this matrix, and cot (a), cot (/3) are the ratios 
of non-diagonal terms to the difference of diagonal ones. 
For example: cot(a) = 2Bs+,s-/{Ks- — Ks+). 

The system of RG differential equations is solved by 
means of an iterative method. 

If cot(a) (orcot(/3)) becomes very large during the flow, 
we stop the flow at some point, introduce the tangent 
of the angles instead of the cotangent, and then resume 
the iteration scheme. In this way, we are able to isolate 
divergences of the prefactors in some of the cosine terms, 
which cause gaps to open and affect observables. 



APPENDIX B: DERIVATION OF THE RG 
EQUATIONS 

1. Flow of the diagonal basis 

To second order in perturbation, one finds the correc- 
tions dKi, dK2, dK^, dK^ to the LL parameters, and 
the non-diagonal terms dBi2, dB^^. These non-diagonal 
terms signal that, after the RG step, Bo is no longer a 
diagonal basis. We then go back to the S-i basis, using 



the transformation 5*"^. In this basis, off-diagonal terms 
have been incremented by small amounts during the RG 
step. For instance 



dB 



1 



s — s-\- 



{dKi - dK2) cos 2a + dBi2 sin 2a (Bl) 



Diagonal terms also undergo infinitesimal variations 
dKs-(s+) = -{dKi+dK2)±-{2dBi2 cos2a+{dKi-dK2) sin 

Similar expressions hold for the charge modes when we 
perform the substitutions s— >c, 1— !-3, 2— +4, a— +/3 
in the equations above. 

The new matrix is diagonalized by the operator S{a + 
da,P + d(3), where the angle da, which account for the 
dB^-^+ and di<r^_(^4.) variations (/i = c, s), indicate a 
rotation of Bo- This idea is summarized in the diagram 
shown in Fig. [3l 

We now determine the renormalization flow of the an- 
gles a and (3. In the spin sector, the diagonalization 
condition is written in terms of and Bs-s+ 



-^{Ks- - Ks+) cos 2a + Bs-s+ sin 2a = (B3) 

One differentiates the above equation in order to relate 
da, and Bs_s+: 



d2a 



sin^ 2a K^- - Ks+ 

\{dKs- — dKs+) cos 2a — dBs-s+ sin 2a 
(sin^ 2a — cos^ 2a) sin 2a 

In the diagonal basis this equation reads 



(B4) 



1 X 2 sin 2a cos 2a 
dcot2a = -— —((dKi-dK2)—^^ 



Ki-K: 



sin 2a — cos^ 2a 

— sin^ 2a + cos^ 2a , , 

dBi2—^^ J—) (B5) 

sm 2a — cos^ 2a 



where the differentials of the LL parameters are known in 
the diagonal basis. They were obtained to second order 
in perturbation, and the dK^, which we use here, were 
given in Sec. IIII Al In the charge sector, we obtain the 
equivalent set of equations with the changes s — > c, 1 — > 
3, 2 — > 4, a ^ /3. The additional expressions for the 
differentials of off-diagonal terms are obtained in a similar 
way, giving, for the case of a generic filling 

dBu^PiQiiigl + al + G^t)- 

KiK2{gla + glc + ale + Gl)) - K^K2h{P,)g,g2 



dB^^ P2Q2{9!c + 9L + 9f\c 



where h{Pi) = {{PiQif + 0.25(p2 _ Ql))-\ 



(B6) 
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2. First order correction to gi, g2 

Setting 511 = gid + gw and 512 = gid - 9id', the Hamil- 
tonian reads 



H = Hll + gi J circos(20s_) cos(20s+)+ 

52 j dr sm{2^s-)sm{2^s+) (B7) 

In order to simplify the RG calculation, we first solve 

this problem in the total/transverse basis where the aver- 
ages over the high energy terms are (0s_(r)^)?i = Ks-dl, 
{<f>s+{rf)h = Ks+dl and {c^s-{r)^s+{r))h = \Adl. One 
may determine the rcnormahzation flow that is produced 
when integrating out the high energy components; for in- 
stance, the renormalization of gi gives 

{gi j drcos{2{(j)s+ + h[(t)s+]))-cos(2{(j)s-+h[(l)s-\)))h = 

j drcos{2{(f>,+ + (f>,_) + {h[(t>s+] + h\(l)s-])))h+ 

+ { J drcos{2{cj),+ - cj),_) + {h[ct>,+] - h[ct>,_])))h (B8) 

We reexponentiate the cosines, use Debye- Waller type 
relations and expand the exponential function in Taylor 

series 



(cos(x + h[x]))h = cos(a;)(^^ exp(2icr/i[x]))/j 

(7=± 

= cos(a;) exp{-{h[xf)h) = (1 - {h[xf) h)cos{x) (B9) 



where ± h[,^,-]))^)h = {K,+ + K,_)dl ± Adl. 

One then finds the usual diagonal term 

fifi(l - {Ks+ + Ks-)dl) j dr{cos{2{(j>s+ + (t>s-)) 

+ cos(2((^s+ - 4>8-)) = 
g,{l - {K,+ + K,_)dl)x 



X j dr cos{2(l)s+) cos{2^s-) (BIO) 



After rcscaling the integration variable dr one gets the 
RG equation for gi. But in the non-diagonal basis there 
is also an additional term 



giAdl j dr(cos(2(^s+ -|- - cos(2((/)s+ - (f>s-)) = 
giAdl j drsin(2(/)s_)sin(20s+) (Bll) 



This links the change of 52 to the coupling constant gi. 
The derivation of the RG equation for the g2 term is 
obtained in a similar fashion, using the identity: 



sin(2(/)s+)sin(2^s_) = ^{cos{(j)s+-(l)s-)-cos{(j)s++(j)s-)) 

(B12) 

Finally, the first-order RG equation for gi is 



^ = 5i(2)(2 - {K,+ + K,-)) + 52(i)A (B13) 
In the diagonal basis, using 

Ks-i2) ^ P'i^l(2) + 2PQA,2 + 0^/^2(1) 

.2 



A ^ PQ{Ki - K2) + Ai2{P' 



The RG equations for the couplings are 







dgi{2) 



= gi(2){K2+K^)+g2(i){{P-Q')A^2+PQ{Ki-K2)) 

(B15) 



3. RG equations for the half filled case 

Using the same method as for the incommensurate case 
we find the following system of equations 
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dKi _ l^jj2t„2 

dl 



dK2 
dl 



+ l{g! + g2) + f{Pi)i9i92))] 

= liQMa + MS)gL + 9l + - K^iPigfa + Phi + M^)Pl9l\\ + QIgI + Mmlglb + Qlgl 

+ \{9l + 9l)-f{Pi){9i92))] 



dKj _ 1 „2r„2 

dl 

dKi 



dl 

dgic 

dl 
dgia 

dl 

dg2c 
dl 
dgi 
dl 
dm 
dl 

dg\\c 
dl 

dg4a 
dl 

dGp 
IT 
dGt 

dl 
dgz\\ 

dl 
dgsa 

dl 

dgzb 
dl 

dgsc 
dl 



^PMc + 9lc + 5yc + 9lc] + l^Ql{9l\\ + gL + ffsfc + 9lc) ' M^) 

■■ \qI\91c + 52c + 9lc + 53c] + \P2{9l\\ + 53a + 536 + 53c) " M^) 

-. gi42 - {PIK2 + PiK^^ + QlKi + QlK^^)] - {gig^c + giag\\c + MS)g3cg3\\) 

= 5la[2 - {PHK2 + K-') + Ql{K, + K^'))] - (3ie5||c + ^0(5)53a53||) 

= 52c[2 - {P\PiK^^ + P^Ki + QlK^^ + QIK2)] - (5lc5l + Jo(<5)53c53a) 
= (?i(2 - {K2 + Ki)) + PiQi{K2 - Ki)g2 - 7(5ic52c + Jo(<5)53653||) 
= 52(2 - {K2 + Ki)) + PiQi{K2 - K,)gi 

= 5||c(2 - iP?K^' + QlK^' + P^K^' + QlK^')) - {giagic + Jo(5)53a53c) 

= 54a(2 - liPhKi + K^') + Ql{K2 + K^'))) 
= - {P^K, + QiK2)) + gliPliKi - i^fi) + QliK2 - K^')) 
= Gt{l - {PlK^' + QlK^')) + gl{P?{-K, + K^') + Q\{-K2 + i^^"')) 
= 5311 (2 - {PIK2 + PlKi + Q\Kr + Q\K^)) - {gigsb + 5ic53c + 5ia53a) 

= 53a(2 - (Pf i^r' + P^Ki + QlK^' + QlKs)) - (5||e53c + 5la53||) 
= 336(2 - (P^Ki + P^Ki + QIK2 + QIK3)) - (gi53|| + 52c<?3c) 

: c/3c(2 - {Pi{K4 + K^^) + QKifa + K^^))) - {g\\^g3a + 52c536 + 5ic53||) 

I 



(B16) 



where 



/(Pi) = (PiQi 



1P2 

4 



PiQ: 



(B17) 



The renormaUzation of the parameter 7 is controlled 
by the same equation as before. The additional flows for 
the velocities of the modes, due to umklapp scattering, 
are all proportional to a Bcsscl term J2(4(^), and hence 
neglected. The general formula describing the flow of the 
diagonal basis remains the same as for the incommensu- 
rate case, but one needs to substitute modified expres- 
sions of the dKy. 



APPENDIX C: ORDER PARAMETER 
OPERATORS IN BOSONIZATION LANGUAGE 



We first write the order parameters in fermionic lan- 
guage. We only consider those order parameters which 

can produce power-law decays of correlations for the var- 
ious locked phase fields combinations. These operators 
are first defined for each site, then expressed in the o/tt 
basis where the Xma coefficients enter their expressions. 

There are two kinds of order parameter operators. The 
first group represent charge density (particle-hole) fluc- 
tuations with a 2 kp wave vector. They correspond to 
the usual CDW, which is the sum of CDW in each band. 
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Up to an unimportant constant factor it gives 

The subscript ^ denotes the band, so, to obtain the or- 
der parameter inside one specific band it is enough to 
take the first or the second term in the above sum. It 
is also possible to define an operator which describes the 
difference of the densities on the two legs 

aa' 

or the operator which describes an orbital antiferromag- 
netic (OAF) fluctuation where currents flow along the 
legs and the rungs of the ladder 

OoAF ■ a^-p,a^aa'a+fj,a' (C3) 

^ aa' 

We can also define (at half filling) "bond" operators, 
which represent density waves located on the bonds, ei- 
ther in phase 

Obdw X! X! exp(fc_F^a;)aL^^(5crcr'a+pcr' (C4) 

out of phase between the two legs of the ladder, 

Ot,bdw X! X! '2xp(fcF^x)a!L^CT'^<T<T'Q:+A'cr' (C5) 

aa' 

or in the diagonal direction: 

Ofdw ~ X! X! exp(fc_F^a;)^ • Q^_^^^b„„'OLJ^^:,a' (C6) 

aa' 

Away from half-filling, on-site and bond operators are de- 
generate, because of translational invariance (the charge 
symmetric mode is masless). 



The second group describes superconducting pairing 
(particle-particle) fluctuations with zero wave vectors. 
As usual there is the s-wave pairing 

Oscs ~ X! X! f^"-/iff'5<T<T'a+pcr' (C7) 

^ aa' 

and the d-wave pairing, which corresponds to a change 
of sign of the order parameter when moving from along 
the legs to along the rungs 

ii aa' 



These phases are given different names in the liter- 
ature. The name orbital antiferromagnet (OAF) was 
used traditionally for the operator defined above, but it 
is also called staggered flux^^ (SF) or d-density wave^ 
phase (DDW). Its bond counterpart is sometimes called 
f-density wavc'^^ (FDW), or diagonal current^^ (DC) 
phase. Similarly, our nCDW and ttBDW orders, are also 
denotedSa CDW and PDW or CDW and SP in Ref. M- 
We have decided to use the notation ttCDW to avoid any 
confusion with the usual CDW which also appears in our 
calculation. 

We can now represent the operators in terms of boson 
fields, using the mapping ([5]). It is important to keep 
the same convention for the signs of the Klein factors as 
that we used to write the Hamiltonian in bosonic form. 
Choosing F — +1 we get 7]a+Va- = This determines 
whether a sin or a cos appears in the formulas below. 
This choice was used in Ref. 2^ and Ref. 46 but the op- 
posite one was used in Ref. [S^. One can easily relate the 
two by shifting the phase fields (fi by an amount tt/2. 

The operators take the form 



OttCdw ~ cos (j)c+ sin 9c- cos 0^+ cos 6s- — sin 0c+ cos dc- sin (/)s+ sin Os- 
OoAF ~ cos 0C+ cos 9c- cos COS 9s- + sin (/)c+ sin Oc- sin (j)s+ sin 9s- 

OnB DW ~ COS 0C+ COS 9c- sin (j)s+ sin 9s- + sin (/I)c+ sin 9c- cos (j)s+ cos 9s- 
O F DW ~ COS (j)c+ sin 9c- sin sin 9s- ~ sin 0c+ cos 9c- cos (j)s+ cos 9s- 
Oscs ^ exp i9c+ cos 9c- sin (j)s+ sin (f>s- — i exp i9c+ sin 9c- cos (j)s+ cos (ps- 
Oscd ^ exp i9c+ cos 9c- cos (l)s+ cos (j>s- — i exp i9c+ sin 9c- sin 0s+ sin (f>s- 



It is also useful to consider these operators in the -Bq/tt are 
basis For example, the SDW operator in the o band 

SDWio) and the CDW operator in the tt band CDWiTi) ^ im ^ a \ ■ 



"^'"^ (CIO) 
Ocdw{7t) ~ exp«((/)c+ - (f>c-) coa{(l)s+ ~ 4)s-) 
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To determine the phases, we need to obtain the expo- 
nents that characterize the spatial decay of the operators' 
correlations . Using a standard procedure to compute the 
correlations with the quadratic Hamiltoniani we find 



VCDW =K2+Ki+Ki + K3 
VOAF = Pi^K2 + QfKi + P*^K^^ 

+ K7^ 



Vscd = K2+K1 



(Cll) 

The exponents of the OAF and ttCDW fluctuations are 
I 



the same, so we need to evaluate logarithmic corrections 
to the powerlaw decay to determine the dominant order- 
ing. 



APPENDIX D: SIMPLIFIED SYSTEM OF RG 
EQUATIONS 

When cot 2a — > and cot 2/3 ^ 00 one gets the follow- 
ing system of first order RG equations for the couplings 



dgic 
dl 


= gic[2 


- iPiK2 + \k^' + 


Ik,')] 


dgia 
dl 




-{pnK2+K^'))] 




dg2c 
dl 


= 52c [2 






dgi 
dl 


= .9i(2- 


- {K2 + K,)) 




dg2 
dl 


= 52(2- 


- {K2 + K,)) 




dg\\c 
dl 


= 5||c(2 


- {PlK-^ + Ik-' 





(Dl) 



r 



The zeroth order approximation to the above system is 
obtained using the fact that K4 {K2) is much smaller 
(larger) than one. 

For the (725*1 phase, the relevance of the important 
coupling needs to be checked. This gives us only one 
differential equation in this case (assuming that close to 
the fixed point [51 1 = jgal = g) 



Taking into account the fact that PiQi < 0, that it 
keeps decreasing during the flow, and that the initial K2 
makes g irrelevant, one finds that a significant decrease 
of K2 would be required in order to make g relevant. 



K,)) + P,Q,{K2-K,)g (D2) 
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